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SUMMARY 

This report describes a numerical method for generating three-dimensional grids for turbomachinery compu- 
tational fluid dynamics codes. The basic method is general and involves the solution of a qusisi-linear elliptic 
partial differential equation via pointwise relaxation with a local relaxation factor. It allows specification 
of the grid point distribution on the boundary surfaces, the grid spacing off the boundary surfaces, and 
the grid orthogonality at the boundary surfaces. It includes adaptive mechanisms to improve smoothness, 
orthogonality, and flow resolution in the grid interior. A geometry preprocessor constructs the grid point dis- 
tributions on the boundary surfaces for general turbomachinery cascades. Representative results are shown 
for a C-grid and an H-grid for a turbine rotor. Two appendices serve as user’s manuals for the basic solver 
and the geometry preprocessor. 


INTRODUCTION 

Three-dimensional computational fluid dynamics codes require computational grids with suitable resolution, 
smoothness, and orthogonality. High grid resolution allows complex flow physics to be modelled near shocks 
and in shear layers. Smoothness of the metric data prevents the flow solution from being dominated by 
truncation error in the metric coefficients. Grid orthogonality at the boundaries simplifies and improves the 
accuracy of any boundary condition involving normal gradients. 

The above qualities are especially difficult to maintain for realistic turbomachinery geometries. Modern three- 
dimensional designs comprise tapered, twisted, leaned, and bowed blade shapes within contoured endwalls 
to tailor secondary flows. Centrifugal compressors and radial turbines involve simultaneous flow turning in 
the meridional and blade-to-blade planes. The periodicity condition within blade rows poses an additional 
unique problem. 

Current grid generation technology is fairly well developed for general two-din)ensional turbomachine cascade 
geometries. Conformal mapping, algebraic interpolation, and partial differential equation methods are all 
used successfully. The general three-dimensional geometry represented by a realistic turbomachine cascade 
has unique requirements for a body-fitted grid that have not been met with current technology. 

This report describes a numerical method for generating three-dimensional grids for turbomachinery com- 
putational fluid dynamics codes. The basic method is general and involves the solution of a quasi-linear 
elliptic partial differential equation via pointwise successive over-relaxation with a local relaxation factor. 
The governing equation contains forcing functions that depend upon the boundary point distribution and the 
boundary surface gradient. The method allows specification of the grid point distribution on the boundary 
surfaces, the grid spacing off the boundary surfaces, and the grid orthogonality at the boundary surfaces. It 
includes adaptive mechanisms to improve smoothness, orthogonality, and flow resolution in the grid interior. 
A geometry preprocessor constructs the grid point distributions on the boundary surfaces for general turbo- 
machinery cascades. It utilizes a two-dimensional version of the basic solver and algebraic interpolation to 
form the boundary distributions for the three-dimensional basic solver. 
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This report includes a description of the coordinate system , a discussion of the mathematical formulation of 
the method, and some representative results. Two appendices serve as the user’s manuals for the geometry 
preprocessor and the basic solver. 


COORDINATE SYSTEM 

The partial differential equations for computational fluid dynamics codes are usually described with reference 
to a generalized coordinate system to simplify the implementation and make it independent of any specific 
geometry. The function of a grid generation system is to generate an ordered distribution of points in 
physical space to align with some body-conforming generalized coordinate system in computational space. 
Physical space is described with reference to a cartesian coordinate system for the grid generation technique 
described in this report. Each boundary surface segment in physical space must coincide with a boundary 
surface segment in computational space. 

Three types of generalized body-conforming coordinate systems are commonly used for turbomachinery 
cascades: C-grid, H-grid, and 0-grid. The grid generation technique described in this report can produce 
C-grids and H-grids. 

Figure 1(a) shows a C-grid about a generic blade shape in the cartesian {xi^X 2 ^x^) coordinate system 
in physical space. The inlet surface is A1-A2-A2'-A1' and the outlet surface is B1-B2-B2'-BP. The hub 
endwall surface is A1'-B1'-B2'-A2' and the shroud endwall surface is A1-B1-B2-A2. The periodic surfaces 
are Al-Bl-Bl'-Al' and A2-B2-B2'-A2^ The blade surface is the wrapped D-E-D-iy-E'-D' surface. Surface 
C-D-D'-C' represents a branch cut from the wrapped blade surface to the outlet surface. Figure 1(b) shows 
the C-grid in the generalized ( 6 )^ 2 )^ 3 ) coordinate system in computational space. 

Figure 2(a) shows an H-grid between two generic blade shapes in the cartesian (a?i,X 2 ,X 3 ) coordinate sys- 
tem in physical space. The inlet surface is A1-A2-A2'-AP and the outlet surface is B1-B2-B2'-B1'. The 
hub endwall is A1'-BP-B2'-A2' and the shroud endwall surface is A1-B1-B2-A2. The periodic surfaces are 
Al-Cl-Cl'-Al', A2-C2-C2'-A2', Dl-Bl-Bl'-DP, and D2-B2-B2'-D2'. The blade surfaces are Cl-Dl-Dl'-Cl' 
and C2-D2-D2'-C2'. Figure 2(b) shows the H-grid in the generalized (^ 1 ,^ 2 , ^ 3 ) coordinate system in com- 
putational space. 


MATHEMATICAL FORMULATION 


The quasi-linear elliptic governing equation is taken from reference 1 as 


3 3 


a d^x dx 


kzzl 


The metric tensor components and 51 ** in equation (1) are defined as 


= a* . aJ 


( 1 ) 

( 2 ) 


where 


a* - dj X dk / y/g L i, k cyclic 
^ = ai • (52 X as) 
dx 


The X in equation (1) represents the position vector in physical space with cartesian (xi, X 2 , X 3 ) components. 
The Pk are the forcing functions specified by the user. They represent one-dimensional stretching in each 
coordinate direction. Values of the forcing functions on the boundaries are determined by specification of 
the boundary point distribution and the boundary surface gradient. Values of the forcing functions in the 
interior are determined by interpolation of the values on the boundaries. 
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Equation (1) can be rewritten using matrix notation as 


A + B P = 0 


( 3 ) 


where 


A = 


/ 9 ^ 

21 

9 

. -31 dH 

\9 d(idU 


«i2 _a!i_ 

-22 2!|, 

9 

„32 d^i 

9 6(ia(3 


-13 _2!I_ 

9 ^^& 


B = 


0 
0 


9 


22 ^ 


0 

0 


p = 



Equation (3) can be solved for P* on the boundaries as 

Po = — Bq^ Ao 


(4) 


where the subscript ‘‘0” indicates values on the boundary. Tangential derivatives for terms on the right hand 
side of equation (4) are determined by applying standard difference formulas to the prescribed boundary 
point distribution on the surface. Normal derivatives are determined by specifying the first normal derivative 
equal to the desired spacing off the boundary and using the approximation 


d^xo _ 2(ii - xo) 2 dxo 
dn? (An)2 An dn 


( 5 ) 


where n indicates the normal direction, the subscript “0” indicates values on the boundary, and the subscript 
“1” indicates values one point away from the boundary. 


Once the boundary values are known, the interior values of P* can be determined using 

3 3 

Pfc(^l,^2,^3) = ^ (1 -f Tib, 1,1 ^Ok.1.3 (1 + «Jb,l,2) Tfc,/,2 

/=:! l-l 


(6) 


where 


a*, 1,1 = C'a*,,., (6 -iUin) 

at, (.2 = C'o*,,,, 

0 k, 1,1 = [( 6 „„ 

0 k,l ,2 = [(^/«„ ~ 

Tfc.I.l = ^fc.l.l/ f Xw 0>i,l,l+^ 0k,l,^ 


r 3 


m=l 

3 


Tfc,/,2 = 13k, 1, 2/ I ^2 4- ^2 I • 


\m=l 


171=1 


The value of Po*.,,! represents the Ar-th component of the Po vector on the minimum boundary surface in 
the /-th direction. The value of Pofc,i,a represents the ib-th component of the Pq vector on the maximum 
boundary surface in the /-th direction. The functions a, and 7 have subscript notation similar to that of 
Pq. The a function represents linear extrapolation from a controlled boundary using a constant factor 
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The 7 function represents the combined effect of the /? functions, which represent power-law factorization 
with constant exponent to control the depth of influence away from a controlled boundary. 

The values of the forcing functions can be modified for improved smoothness by using 

Pk = SP, (7) 


where 

0=1- tanh(C<?j(l — 

The constants C$^ and Cb^ define the rate and order of the adaptation. The variable is a measure of the 
shear of the grid, with J representing the Jacobian and liyhjh representing the grid cell lengths in each 
direction. 

A measure of the local orthogonality of the grid can be defined as 


(j> — jj. 


( 8 ) 


The constant exponent defines the order of the adaptation. A one-dimensional variational form of 
equation (8) can then be written as 

= (9) 


If a flow variable gradient E is computed by the flow solver such that £■ > 0, a measure of the local flow 
resolution can be defined as 

V^ = J(l + ^). (10) 

A one-dimensional variational form of equation (10) can then be written as 


i>'k = 


h{J) 

J 


6k{E) 

1 + E' 


( 11 ) 


The values of the forcing functions can be modified for improved local orthogonality and flow resolution 
using 

P^ = (l + At)P* F,.P,>0 

P;, = il-\k)Pk Fk-Pk<0 ^ ’ 

where 

Xk = Cx, Unh(Cx,\Fkf^^) 

Fk = (W^<l>'k + W^^'k)/(W^ + W^) 

with the constants C \^ , and C\^ determining the range, rate, and order of the adaptation. The variable 
F represents a weighted combination of the skewness and flow error variations, with and W^p as the 
respective weighting constants. 


These adaptive adjustments to the forcing functions can be employed in an accumulative or a non- 
accumulative manner. Under the accumulative method, the adjustments are lagged, but the forcing functions 
always satisfy the physical constraints. Under the non-accumulative method, the adjustments are immediate, 
but they are limited by the range and rate constants. 


The following sequence is iterated until convergence is attained: 

1. Solve the governing equation with current forcing functions using pointwise relaxation with a local 
relaxation factor for stability. 

2. Evaluate the second normal derivatives on the boundary surfaces. 

3. Evaluate the forcing functions on the boundary surfaces. 
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4. Evaluate the forcing functions in the interior by interpolation. 

5. Adjust the forcing functions for adaptive smoothness, orthogonality, and flow resolution. 

RESULTS 

Figure 3 shows a three-dimensional C-grid for a turbine rotor. The grid comprises 101 points in the 
direction around the blade, 16 points in the ^2 direction away froni the blade, and 15 points in the ^3 
direction from hub to tip. The grid was truncated in the direction and thinned in the ^2 direction to 
improve the appearance of the figure. The spacing at the solid surfaces was specified as 0.1 relative to the 
uniform unit spacing. Default values were used for all other variables. This grid required approximately 
1200 seconds of Cray X-MP computing time for 600 iterations. 

Figure 4 shows a three-dimensional H-grid for the same turbine rotor as Figure 3. The grid comprises 61 
points in the streamwise direction, 31 points in the pitchwise (2 direction, and 15 points in the radial (3 
direction. The grid was thinned in the ^2 direction to improve the appearance of the figure. The spacing 
at the solid surfaces was specified as 0.1 relative to the uniform unit spacing. Default values were used for 
all other variables. This grid required approximately 1400 seconds of Cray X-MP computing time for 600 
iterations. 


CONCLUDING REMARKS 

A general numerical method for generating three-dimensional grids was developed and implemented along 
with a geometry preprocessor for turbomachinery cascades. Results were shown for a C-grid and an H- 
grid for.^a turbine rotor. The method includes an adaptive mechanism for improved flow resolution when 
coupled with a flow solver, but this was not demonstrated. Since the basic method is completely general, 
additional preprocessors for other physical geometries could be developed to extend the application of this 
grid generation method. 
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APPENDIX A 


USER’S MANUAL FOR BASIC SOLVER 

The basic solver is coded in FORTRAN IV as program GRID3D with 34 subprograms. The nonstandard 
PARAMETER statement is used to facilitate dimensioning of arrays, while the nonstandard NAMELIST 
feature is used for quasi-free-format input. Figure 5 shows a flow chart of the program. 

Input to the basic solver comprises three namelists read on unit 5 and the initial grid binary file read on 
unit 10. The type (REAL or INTEGER) of each variable follows standard FORTRAN convention. The I, 
J, and K indices correspond to the and ^3 directions, respectively. 

Namelist SYSTEM 


ISYS 

= 1 

Solve governing equation as a Laplacian system without forcing functions (default). 
This option produces a grid that tends toward uniform spacing away from curved 
boundaries, but produces smaller spacing near convex boundaries and larger spacing 
near concave boundaries. 


= 2 

Solve governing equation as a Poission system with forcing functions. This option 
allows specificaton of spacing and orthogonality at the boundaries. 

INTERP 

= 1 

Run three-dimensional interpolation as initial guess (default). This option is used 
when the geometry preprocessor is employed to generate the initial grid binary 
file without interior values. 


= 0 

Bypass three-dimensional interpolation. This option is used when a complete 
initial grid binary file is imported from another grid generator. 

ICHECK 

= 1 

Check boundary structures for grid folding and write report on unit 20. 


= 0 

Bypass boundary structure check (default). 

ICONT 

= 0 

Terminate operation if boundary structure grid folding is detected. 


= 1 

Continue operation (default). 

ITMAX 


Maximum number of overall iterations (default = 100). 

LREFI 

= 1 

Establish a reference length scale for spacing off the 1 = 1 and I = IMAX surfaces 
based upon the unit length in the J-direction for each K-layer at the 1=1 surface. 
The unit length is defined as (total length)/(number of points - 1). 


= 2 

Establish a reference length scale for spacing off the 1=1 and I = IMAX surfaces 
based upon the local unit length in the J-direction. The unit length is defined as 
(total length)/(number of points - 1 ). 


= 3 

Establish a reference length scale for spacing off the 1=1 and I = IMAX surfaces 
based upon the value of FACIl from the initial grid binary file or the 
three-dimensional interpolation (default). 

LREFJ 


Same as above for spacing off the J = 1 and J = JMAX surfaces. 

LREFK 


Same as above for spacing off the K = 1 and K = KM AX surfaces. 

RLXSOR 


Relaxation factor for the SOR solver. Suitable values range from 0.0 to 2.0 
(default = 1 . 0 ). 

RLXBDE 


Relaxation factor for extrapolated second derivatives at the boundary surfaces 
(default = 0.15). 

RLXADP 


Relaxation factor for grid adaptation (default = 0.15). 

EPSFLD 


Folding grid point criterion (default = 0.0). 

EPSBAD 


Sheared grid point criterion (default = sine 10®). 

EPSCNV 


Convergence criterion (default = 0.0001). 

Namelist PBASIC 


IlCONT 

= 1 

Control spacing and orthogonality off the 1=1 surface. 


= 0 

Do not control spacing and orthogonality off the 1=1 surface (default). 
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IMCONT 

JICONT 

JMCONT 

KICONT 

KMCONT 

DELTIl 

DELTIM 

DELTJl 

DELTJM 

DELTKl 

DELTKM 


Same as above for I = IMAX surface. 

Same as above for J =r 1 surface. 

Same as above for J ~ JMAX surface. 

Same as above for K = 1 surface. 

Same as above for K = KMAX surface. 

Desired spacing off the 1=1 surface. If DELTIl = 0.0, the value will be set equal 
to the value of EACH from the initial grid binary file (default). 

Same as above for the desired spacing off the I = IMAX surface. 

Same as above for the desired spacing off the J = 1 surface. 

Same as above for the desired spacing off the J = JMAX surface. 

Same as above for the desired spacing off the K = 1 surface. 

Same as above for the desired spacing off the K = KMAX surface. 


The following parameters control the depth of grid clustering away from the controlled surface. The 
default value of 0.0 produces a quasi-linear rate of spacing increase away from the boundary. Positive 
values produce a nonlinear rate with continuously increasing spacing in the interior. Negative values 
produce a nonlinear rate with near-uniform spacing in the interior. 


APIl 


APIM 


APJl 


APJM 


APKl 


APKM 

^01.3, 3 

AQll 


AQIM 

^«3,1,3 

AQJl 

^03,3,1 

AQJM 


AQKl 

^03,3,1 

AQKM 

^Of3,3,3 

ARIl 


ARIM 

^03,1,2 

ARJl 

^Cr3.2,l 

ARJM 

^CT3,2,2 

ARKl 

^ora.3,1 

ARKM 

^Cr3,3,2 


for -direction extrapolation 
for -direction extrapolation 
for ^ 2 'direction extrapolation 
for ^ 2 'direction extrapolation 
for ^ 3 -direction extrapolation 
for ^ 3 -direction extrapolation 
for -direction extrapolation 
for ^i-direction extrapolation 
for ^ 2 'direction extrapolation 
for ^ 2 -direction extrapolation 
for ^ 3 -direction extrapolation 
for ^ 3 -direction extrapolation 
for ^i-direction extrapolation 
for -direction extrapolation 
for ^ 2 -direction extrapolation 
for ^ 2 -direction extrapolation 
for ^ 3 -direction extrapolation 
for ^ 3 -direction extrapolation 


of Pi from 1=1 surface, 
of Pi from I = IMAX surface, 
of Pi from J = 1 surface, 
of Pi from J = JMAX surface, 
of Pi from K = 1 surface, 
of Pi from K = KMAX surface, 
of P 2 from 1=1 surface, 
of P 2 from I = IMAX surface, 
of P 2 from J = 1 surface, 
of P 2 from J = JMAX surface, 
of P 2 from K = 1 surface, 
of P 2 from K = KMAX surface, 
of P 3 from 1=1 surface, 
of P 3 from I = IMAX surface, 
of P 3 from J = 1 surface, 
of P 3 from J = JMAX surface, 
of P 3 from K = 1 surface, 
of P 3 from K = KMAX surface. 


The following parameters control the depth of grid orthogonality away from the controlled boundary. 
Larger positive values produce incresised depth (default = 3.0). 


BIl 

J for ^ 1 -direction 

BIM 

3 for ^i-direction 

BJl 

, for ^ 2 -direction 

BJM 

3 for ^ 2 -direction 

BKl 

C/J 3 J for ^ 3 -direction 

BKM 

3 for ^ 3 -direction 


power-law factorization of Pi,P2, P3 
power-law factorization of Pi, P2, P3 
power-law factorization of Pi,P2, P3 
power-law factorization of Pi, P2, P3 
power-law factorization of Pi,P2, P3 
power-law factorization of Pi, P2, P3 


from 1=1 surface 
from I = IMAX surface 
from J = 1 surface 
from J = JMAX surface 
from K = 1 surface 
from K = KMAX surface 


Namelist PADAPT 


SMRATE 

SMORDR 

SMBDWT 


Rate of the penalty function for smoothness adaptation (default = 3.0). 

Larger values produce smoother grids tending toward the Laplacian solution. 

Power of the penalty function for smoothness adaptation (default = 1.0). 

Larger values produce more smoothing near discontinuous areas. 

Parameter providing boundary protection for smoothness adaptation (default = 0.0). 
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Larger values allow deeper penetration of the specified boundary conditions into the 
interior at the expense of smoothness adaptation. The maximum value of unity 
produces a linear penetration. 

OFRNGE Range of the penalty function for combined orthogonality and resolution adaptation, 

when lACCUM = 0 (default = 1.0). Larger values produce more combined adaptation. 
Lag of the combined orthogonality and resolution adaptation when lACCUM = 1. 

A zero value produces no lag, while a unity value produces full lag. 

OFRATE Rate of the penalty function for combined orthogonality and resolution adaptation. 

Suggested values are OFRATE = 5.0 when lACCUM = 0 and OFRATE = 50.0 
when lACCUM = 1. Larger values produce more combined adaptation. 

OFORDR Power of the penalty function for combined orthogonality and resolution (default = 1.0). 

Larger values produce more adaptation near highly sheared and high gradient areas. 

OFBDWT Parameter providing boundary protection for combined orthogonality and resolution 

adaptation (default = 0.0). Larger values allow deeper penetration of the specified 
boundary conditions into the interior at the expense of combined adaptation. The 
maximum value of unity produces a linear penetration. 

OTORDR Power of the skewness function for orthogonality adaptation (default = 2.0). 

Larger values produce more adaptation near highly sheared areas. 

WTOT Weight factor for relative effect of orthogonality adaptation (default = 1.0). 

WTFW Weight factor for relative effect of resolution adaptation (default = 0.0). 

lACCUM = 1 Perform accumulative adaptation. The adjustments to the forcing functions 

are lagged, but the forcing functions always satisfy the physical constraints (default). 

= 0 Perform non accumulative adaptation. The adjustments to the forcing functions 
are immediate and are limited only by the range and rate constants. 

The initial grid binary file is read from FORTRAN unit 10 with the following code: 

READ (10) IMAX»JNAX,KNAX 

DO 100 K=1,KMAX 

READ (10) ((X(I,J,K),I=1,IMAX).J=1,JMAX). 

* ((Y(I,J,K),I=1,IMAX),J=1,JMAX), 

* ((Z(I,J,K),I=1,IKAX),J=1,JMAX) 

100 COHTINUE 

READ (10) IGTYPE 

READ (10) ISING1,ISING2,ISING3,ISIKG4 

READ (10) JSING1,JSING2,JSING3,JSIKG4 

READ (10) FACI1,FACIM 

READ (10) FACJ1,FACJM 

READ (10) FACKl.FACKM 

The variables in the binary file are defined as follows: 

Number of grid points in the direction. 

Number of grid points in the ^2 direction. 

Number of grid points in the ^3 direction. 

Physical Cartesian coordinate in the xi direction. 

Physical Cartesian coordinate in the X 2 direction. 

Physical Cartesian coordinate in the X 3 direction. 

1 C-grid. 

2 H-grid. 

I-index of first line singularity on the J = 1 surface. 

I-index of second line singularity on the J = 1 surface. 

I-index of first line singularity on the J = JMAX surface. 

I-index of second line singularity on the J = JMAX surface. 


IMAX 

JMAX 

KMAX 

X 

Y 

Z 

IGTYPE = 

ISINGl 

ISING2 

ISING3 

ISING4 
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JSINGl J-index of first line singularity on the I = 1 surface. 

JSING2 J'index of second line singularity on the 1=1 surface. 

JSING3 J-index of first line singularity on the J = JMAX surface. 

JSING4 J-index of second line singularity on the J = JMAX surface. 

The following values are only used if DELTIl, DELTIM, DELTJl, DELTJM, DELTKl, and DELTKM 
in namelist PBASIC are defaulted. 

EACH Desired spacing off the 1=1 surface. 

FACIM Desired spacing off the I = IMAX surface. 

FACJ 1 Desired spacing off the J = 1 surface. 

FACJM Desired spacing off the J = JMAX surface. 

FACKl Desired spacing off the K = 1 surface. 

FACKM Desired spacing off the K = KMAX surface. 

Output from the basic solver consists of five files. The system message file is written on unit 6. The locations 
of any boundary and interior folding points are written on units 20 and 30, respectively. The final grid file 
is written on unit 40 with the same format as the initial grid file. A plot file for graphics post<processing is 
written on unit 50. 

Folding points are identified where the normalized Jacobian is negative. Sheared points are identified where 
the normalized Jacobian is less than EPSBAD. An error index is computed every ten iterations and is defined 
as the relative movement of a point normalized by the diagonal length of the grid. Convergence is attained 
when the absolute maximum error index among all points is less than EPSCNV with no folding points. 
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APPENDIX B 


USER’S MANUAL FOR GEOMETRY PREPROCESSOR 

The preprocessor is coded in FORTRAN IV as program BLADE with 39 subprograms, A two-dimensional 
version of the basic solver is included as subroutine GRID2D to define nodal point distributions on the hub 
and shroud boundary surfaces. The nonstandard PARAMETER statement is used to facilitate dimensioning 
of arrays, while the nonstandard NAMELIST feature is used for quasi-free-format input. 

Input to the preprocessor comprises three namelists read on unit 5 for GRID2D, one namelist read on unit 7 
for BLADE, and formatted blade geometry data on unit 8. The type (REAL or INTEGER) of each variable 
follows standard FORTRAN convention. The I, J, and K indices correspond to the ^2> and ^3 directions, 
respectively. 

Namelist SYSTEM (Unit 5) 

ISYS == 1 Solve governing equation as a Laplacian system without forcing functions (default). 

This option produces a grid that tends toward uniform spacing away from curved 
boundaries, but produces smaller spacing near convex boundaries and larger spacing 
near concave boundaries. 

= 2 Solve governing equation as a Poission system with forcing functions. This option 
allows specificaton of spacing and orthogonality at the boundaries, 

ICHECK = 1 Check boundary structures for grid folding and write report on unit 20. 

= 0 Bypass boundary structure check (default). 

ICONT = 0 Terminate operation if boundary structure grid folding is detected. 

= 1 Continue operation (default). 

ITMAX Maximum number of overall iterations (default = 100). 

RLXSOR Relaxation factor for the SOR solver. Suitable values range from 0.0 to 2,0 

(default = 1.0). 

RLXBDE Relaxation factor for extrapolated second derivatives at the boundary surfaces 

(default = 0.15). 

RLXADP Relaxation factor for grid adaptation (default = 0.15). 

EPSFLD Folding grid point criterion (default = 0.0). 

EPSBAD Sheared grid point criterion (default = sine 10®). 

EPSCNV Convergence criterion (default = 0.0001). 

Namelist PBASIC (Unit 5) 

IlCONT = 1 Control spacing and orthogonality off the I =1 surface. 

= 0 Do not control spacing and orthogonality off the 1=1 surface (default), 

IMCONT Same as above for I = IMAX surface. 

J ICONT Same as above for J = 1 surface. 

JMCONT Same as above for J = JMAX surface. 

DELTIl Desired spacing off the 1 = 1 surface. If DELTIl = 0.0, the value will be set equal 

to the value of FACIl from the initial grid binary file (default). 

DELTIM Same as above for the desired spacing off the I = IMAX surface. 

DELTJ 1 Same as above for the desired spacing off the J = 1 surface. 

DELTJM Same as above for the desired spacing off the J = JMAX surface. 

The following parameters control the depth of grid clustering away from the controlled surface. The 
default value of 0.0 produces a quasi-linear rate of spacing increase away from the boundary. Positive 
values produce a nonlinear rate with continuously increasing spacing in the interior. Negative values 
produce a nonlinear rate with near-uniform spacing in the interior. 
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APIl 

^^1,1, 1 

APIM 


APJl 


APJM 

^«1,3,3 

AQIl 


AQIM 


AQJl 

^03,2,1 

AQJM 

^«3,3,3 


for ^ 1 -direction extrapolation of Pi from 1=1 surface, 
for 4i-direction extrapolation of Pi from I = IMAX surface, 
for ^ 2 -direction extrapolation of Pi from J = 1 surface, 
for ^ 2 -direction extrapolation of Pi from J = JMAX surface, 
for ^ 1 -direction extrapolation of P 2 from 1=1 surface, 
for ^i-direction extrapolation of P 2 from I = IMAX surface, 
for ^ 2 “direction extrapolation of P 2 from J = 1 surface, 
for ^ 2 -direction extrapolation of P 2 from J = JMAX surface. 


The following parameters control the depth of grid orthogonality away from the controlled boundary. 
Larger positive values produce increased depth (default = 3.0). 

BIl ^i-direction power-law factorization of Pi and P 2 from 1 = 1 surface 

BIM Cpi ,2 for ^i-direction power-law factorization of Pi and P 2 from I = IMAX surface 

BJl for ^ 2 -direction power-law factorization of Pi and P 2 from J = 1 surface 

BJM for ^ 2 -directibn power-law factorization of Pi and P 2 from J = JMAX surface 


Namelist PADAPT (Unit 5) 


SMRATE 

SMORDR 

SMBDWT 

OFRNGE 

OFRATE 

OFORDR 

OFBDWT 

OTORDR 

WTOT 

WTFW 

lACCUM = 1 


Rate of the penalty function for smoothness adaptation (default = 3.0). 

Larger values produce smoother grids tending toward the Laplacian solution. 

Power of the penalty function for smoothness adaptation (default = 1.0). 

Larger values produce more smoothing near discontinuous areas. 

Parameter providing boundary protection for smoothness adaptation (default = 0.0). 
Larger values allow deeper penetration of the specified boundary conditions into the 
interior at the expense of smoothness adaptation. The maximum value of unity 
produces a linear penetration. 

Range of the penalty function for combined orthogonality and resolution adaptation, 
when lACCUM = 0 (default = 1.0). Larger values produce more combined adaptation. 
Lag of the combined orthogonality and resolution adaptation when lACCUM = 1. 

A zero value produces no lag, while a unity value produces full lag. 

Rate of the pensdty function for combined orthogonality and resolution adaptation. 
Suggested values are OFRATE = 5.0 when lACCUM = 0 and OFRATE = 50.0 
when lACCUM = 1. Larger values produce more combined adaptation. 

Power of the penalty function for combined orthogonality and resolution (default = 1.0). 
Larger values produce more adaptation near highly sheared and high gradient areas. 
Parameter providing boundary protection for combined orthogonality and resolution 
adaptation (default = 0.0). Larger values allow deeper penetration of the specified 
boundary conditions into the interior at the expense of combined adaptation. The 
maximum value of unity produces a linear penetration. 

Power of the skewness function for orthogonality adaptation (default = 2.0). 

Larger values produce more adaptation near highly sheared areas. 

Weight factor for relative effect of orthogonality adaptation (default = 1.0). 

Weight factor for relative effect of resolution adaptation (default = 0.0). 

Perform accumulative adaptation. The adjustments to the forcing functions 
are lagged, but the forcing functions always satisfy the physical constraints (default). 
Perform nonaccumulative adaptation. The adjustments to the forcing functions 
are immediate and are limited only by the range and rate constants. 


Namelist BLDATA (Unit 7) 


IGTYPE = 1 Generate C-grid. 

= 2 Generate H-grid. 

DELTAT Periodic pitch angle in radians. 
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KCUT 

KMAX 

FACKl 

FACK2 

PWK 

NBLD 

FACBl 

FACB2 

PWBD 

NLEAD 

FLCUTl 


FLCUT2 

NTAIL 

FTCUTl 


FTCUT2 

PWCUT 

NIO 

FACIOl 

FACI02 


Number of blade geometry cuts in the K-direction supplied on unit 8. 

Desired number of grid points in the K-direction. 

Desired spacing off the K = 1 surface. 

Desired spacing off the K = KMAX surface. 

Power of stretching function for the K-direction (default = 1.0). Larger values 
produce tighter clustering near the boundaries. 

Desired number of points on upper and lower blade surfaces. 

Desired spacing on the blade at leading edge. 

Desired spacing on the blade at trailing edge. 

Power of stretching function for the blade point distribution (default = 1.0). 

Larger values produce tighter clustering at the boundaries. 

Desired number of points on the leading edge branch cut (H-grid only). 

Desired spacing on leading edge branch cut at leading edge (H-grid only). 

If FLCUTl > 100.0, the spacing on the branch cut at the leading edge will 
match the spacing on the blade at the leading edge. 

Desired spacing on leading edge branch cut at inlet (H-grid only). 

Desired number of points on trailing edge branch cut. 

Desired spacing on trailing edge branch cut at trailing edge. 

If FTCUTl > 100.0, the spacing on the branch cut at the trailing edge will 
match the spacing on the blade at the trailing edge. 

Desired spacing on trailing edge branch cut at outlet. 

Power of stretching function for the branch cut point distributions (default = 1.0). 
Larger values produce tighter clustering at the boundaries. 

Desired number of points in the J-direction at the outlet from trailing edge branch 
cut to the periodic lines for C-grid or one-half the desired number of points in the 
J-direction at the inlet and outlet for H-grid. 

Desired spacing in the J-direction off the branch cut at the outlet for C-grid 
or desired spacing in the J-direction off the periodic surfaces at the inlet and outlet 
for H-grid. 

Desired spacing in the J-direction off the periodic surfaces at the outlet for C-grid 
or desired spacing in the J-direction at the center of the inlet and outlet for H-grid. 


The following parameters apply only to the C-grid: 


ISL 

ISU 

FPDA 

FPDB 

PWIN 

FPDC 

FPDD 

PWPRD 


Desired number of points on lower half of inlet section. 

Desired number of points on upper half of inlet section. 

Desired spacing at center of inlet section. 

Desired spacing at edges of inlet section. 

Power of stretching function for inlet section point distribution (default = 1.0). 
Larger values produce tighter clustering at the boundaries. 

Desired spacing on the periodic surfaces at the inlet. 

Desired spacing on the periodic surfaces at the outlet. 

Power of stretching function for point distribution on the periodic surfaces 
(default = 1.0). Larger values produce tighter clustering at the boundaries. 


The formatted blade geometry data is read from FORTRAN unit 8 with the following code: 


DO 680 K=1,KCUT 
READ (8,1010) HPATH 


READ (8,1020) 
READ (8,1020) 


(ZPATH(I), 1=1, HPATH) 
(RPATH(I), 1=1, HPATH) 
READ (8,1010) IBLDP 
READ (8,1020) (BLDPZC(I) ,1=1 ,IBLDP) 
READ (8,1020) (BLDPTC(I) , 1=1 ,IBLDP) 


12 


READ (8.1020) 
READ (8.1010) 
READ (8.1020) 
READ (8.1020) 
READ (8.1020) 
READ (8.1020) 


(BLDPRC(I).I=1.IBLDP) 

IBLDS 

(BLDSZC(I).I=1.IBLDS) 
(BLDSTC(I). 1=1. IBLDS) 
(BLDSRC(I). 1=1. IBLDS) 
SLOPL.SLOPT 


680 CONTINUE 


1010 FORMAT (113) 

1020 FORMAT (8F10.0) 

The variables in the formatted file are defined as follows: 


NPATH 

ZPATH 

RPATH 

IBLDP 

BLDPZC 

BLDPTC 

BLDPRC 

IBLDS 

BLDSZC 

BLDSTC 

BLDSRC 

SLOPL 

SLOPT 


Total number of data points defining the surface cut. 

Axial coordinates of data points on the surface cut. 

Radial coordinates of data points on the surface cut. 

Total number of data points defining the pressure or lower surface of blade. 

Axial coordinates of data points on pressure surface. 

Circumferential coordinates of data points on pressure surface (in radians). 
Radial coordinates of data points on pressure surface. 

Total number of data points defining the suction or upper surface of blade. 

Axial coordinates of data points on suction surface. 

Circumferential coordinates of data points on suction surface (in radians). 

Radial coordinates of data points on suction surface. 

Tangent of inflow angle or leading edge mean camber angle. If SLOPL > 100.0, 
SLOPL will be reset to match the computed leading edge mean camber angle. 
Tangent of outflow angle or trailing edge mean camber angle. If SLOPT > 100.0, 
SLOPT will be reset to match the computed trailing edge mean camber angle. 


Output from the preprocessor consists of eight files. The system message file is written on unit 6. The initial 
grid file is written on unit 10 for use by the basic solver. The locations of any boundary folding points on 
the hub and shroud surfaces are written unit 20. The location of any interior folding points on the hub and 
shroud surfaces are written on unit 30. The computed nodal point distributions for the hub and shroud 
surfaces are written on units 40 and 45. Plot files for graphics post-processing are written on units 50 and 55 
for the hub and shroud surfaces. 


Folding points are identified where the normalized Jacobian is negative. Sheared points are identified where 
the normalized Jacobian is less than EPSBAD. An error index is computed every ten iterations and is defined 
as the relative movement of a point normalized by the diagonal length of the grid. Convergence is attained 
when the absolute mziximum error index among all points is less than EPSCNV with no folding points. 

The initial grid binary file is written on FORTRAN unit 10 with the following code: 

WRITE (10) IMAX.JMAX.KMAX 
DO 100 K=1.KMAX 

WRITE (10) ((X(I.J.K).I=1.IMAX).J=1.JMAX). 

♦ ((Y(I.J.K).I=1.IMAX).J=1.JMAX). 

♦ ((Z(I.J.K).I=1.IMAX).J=1.JMAX) 

100 CONTINUE 


13 



WRITE (10) 
WRITE (10) 
WRITE (10) 
WRITE (10) 
WRITE (10) 
WRITE (10) 


IGTYPE 

ISIIGl , ISIIG2 , ISIIG3 , ISIIG4 

JSIIGl , JSIIG2 , JSIIG3 , JSIIG4 

FACIl.FACIM 

FACJ1,FACJM 

FACKl.FACKN 


The variables in the binary file are defined as follows: 


IMAX 

JMAX 

KMAX 

X 

Y 

Z 

IGTYPE = 1 
= 2 

ISINGl 

ISING2 

ISING3 

ISING4 

JSINGl 

JSING2 

JSING3 

JSING4 

EACH 

FACIM 

FACJl 

FACJM 

FACKl 

FACKM 


Number of grid points in the direction. 

Number of grid points in the ^2 direction. 

Number of grid points in the ^3 direction. 

Physical Cartesian coordinate in the X\ direction. 

Physical Cartesian coordinate in the X 2 direction. 

Physical Cartesian coordinate in the X 3 direction. 

C-grid. 

H-grid. 

I-index of first line singularity on the J = 1 surface. 

I-index of second line singularity on the J = 1 surface. 
I-index of first line singularity on the J = JMAX surface. 
I"index of second line singularity on the J = JMAX surface. 
J-index of first line singularity on the 1 = 1 surface, 

J-index of second line singularity on the 1=1 surface. 
J"index of first line singularity on the J = JMAX surface. 
J-index of second line singularity on the J = JMAX surface. 
Desired spacing off the 1=1 surface. 

Desired spacing off the I = IMAX surface. 

Desired spacing off the J = 1 surface. 

Desired spacing off the J = JMAX surface. 

Desired spacing off the K = 1 surface. 

Desired spacing off the K = KMAX surface. 
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(B) C-GRID IN COMPUTATIONAL SPACE. 
FIGURE 1. - C-6RID ABOUT GENERIC BLADE. 
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(A) H-GRID IN PHYSICAL SPACE. 
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(B) H-GRlD IN COMPUTATIONAL SPACE. 
FIGURE 2. - H-GRID BETHEEN GENERIC BLADES. 
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FIGURE 3. - THREE-DIMENSIONAL C-GRID FOR TURBINE ROTOR. 


FIGURE <t. - THREE-DIMENSIONAL H-GRID FOR TURBINE ROTOR. 
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f\Gm S. - FLOW CHART FOR PROGRM 6RID 3D. 



FIGURE 5. - CONTINUED. 





FIGURE S. - CONTINUED. 
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FIGURE 5. - CONTINUED. 




FIGURE 5. - CONCLUDED. 
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